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Abstract 

Starting from the ++ glueball mass and wave function computed from lattice QCD, we compute 
the local potential between two constituent gluons. Since the properties of constituent gluons 
are still a matter of research, we allow for them to be either massless, or massive with a mass 

CN ■ 

around 0.7 GeV. Both pictures are actually used in the literature. When the gluons are massless, 
the corresponding local potential is shown to be compatible with a Cornell form, that is a linear 
confinement plus a short-range Coulomb part, with standard values for the flux tube energy density 
and for the strong coupling constant. When the gluons are massive, the confining potential is a 



saturating one, commonly used to simulate string-breaking effects. These results fill a gap between 
lattice QCD and phenomenological models: The picture of the scalar glueball as a bound state 



of two constituent gluons interacting via a phenomenological potential is shown to emerge from 

pure gauge lattice QCD computations. Moreover, we show that the allowed potential shape is 

■ , constrained by the mass of the constituent gluons. 
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I. INTRODUCTION 



The study of glueballs currently deserves much interest from a theoretical point of view, 
either within the framework of lattice QCD or effective models. Their experimental detection 
is also an active field of research-for recent reviews, see Refs. jl| . Among the various effective 
approaches which have been proposed (bag model [2J, QCD in Coulomb gauge [3(, . . . ), many 
studies have been devoted to potential models of glueballs. In this framework, glueballs are 
seen as bound states of two or more constituent gluons interacting via a phenomenological 
potential. Early works on this subject are quoted in Refs. [4j, [5], and more recent studies 
can be found in Refs. 

Actually, the relevance of using a potential model to describe a glueball is still contro- 
versial. Assuming that a potential model can be used, two basic questions then appear: 
What is the mass of a constituent gluon, and what is the potential? The answer to these 
questions differ from one approach to another. Let us begin by the problem of the mass. On 
the one hand, it is argued in some works that a gluon is a massless particle, which acquires 
a dynamical mass given by a = ( \ffP) because of the confining interaction. Relativistic, 

n n Q 

spin-dependent, corrections, are then developed in powers of l//x 2 [4j, la, u\ ■ In this picture, 
the constituent gluon is a posteriori massive, because it is confined into a glueball, and 
this mass is state-dependent. It is worth mentioning that, more generally, both quarks and 
gluons can acquire a constituent mass from renormalization theory. This constituent mass 
runs with the momentum: One can look in particular at the Coulomb gauge approach of 
Ref. Q, where it is shown that massless gluons acquire a running mass which is about 0.7 
GeV at zero momentum. The same approach can also be applied to compute the constituent 



mass of light quarks llj. On the other hand, it is often assumed in different studies that 
a constituent gluon has to be a priori considered as massive 0, [ijj]. The idea underlying 
these last approaches is roughly that the nonperturbative effects of QCD causes a mass 
term to appear in the gluon propagator. Consequently, the gluons should be seen as massive 
particles, with a fixed mass which is typically assumed to be around 0.5 ± 0.2 GeV j^, \\\ . 
The relativistic corrections are then expanded as usual in powers of Interestingly, the 

typical value of m g in the second approach is compatible with the dynamical mass /i for the 
ground state in the first approach, and with the constituent gluon mass at zero momentum 
of the Coulomb gauge model of Ref. [lol |. 
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We turn now our attention to the potentials appearing in the various existing models 
of QCD. In a two-body system, the best-known phenomenological potential is the Cornell 
one, which is roughly of the form ar — n/r, r being the separation between the confined 
particles. It is worth mentioning that ar is the energy of a straight string of energy density 
a, also called the flux tube, linking the quark to the antiquark and encoding the confining 
interaction. The Coulomb part is the lowest order contribution of the one gluon exchange 
processes. The Cornell potential arises from QCD in the case of a quark-antiquark bound 
state, as it can be shown by the Wilson loop technique [14]]. Lattice QCD computations 
of the energy between a static quark-antiquark pair also support this potential jl5|, p. 42]. 
Furthermore, background perturbation theory tells that the potential between two massless 
constituent gluons should also be of the Cornell form [l6j|. Bound states of gluons with the 
Cornell potential have been investigated for example in Refs. Q, If], 3, Q]- But, as a linearly 
rising potential neglects string breaking effects, which have been observed in lattice QCD 
between static quarks, another confining potential is also often used, that is a saturating 
one of the form 2m g (l — e~ r ^ Vc ). As for the Cornell one, models built on such a potential 
have been applied to usual hadrons [3], but also to glueballs j^, Q. Let us note that the 
short-range part corresponding to massive constituent gluons is not a Coulomb term, but is 
proportional to the Yukawa potential e~ m s r jr. 

/,From this discussion, we can conclude that the best way of dealing with constituent 
gluons is still controversial. Consequently, it is of particular interest to try to obtain relevant 
informations from more fundamental approaches such as lattice QCD. If the mass spectrum 
of pure gauge QCD-the glueball spectrum-is now accurately computed in lattice QCD |l9| . 
the potential energy between two constituent gluons has been much less studied than the 
quark-antiquark one. Up to know, the only method to obtain this energy with lattice QCD 
is to compute the energy between two static sources, these sources being in the adjoint 



representation of SU(3) [15l . p. 69]. The Cornell shape is then favored, as in the quark- 
antiquark case. However, nowadays, both the masses and wave functions of glueballs can 
be computed by lattice calculations 20, 2lj]. We propose in this paper a new method for 



extracting the potential between two constituent gluons from these lattice QCD data. Such 
a method has the conceptual advantage of dealing with "physical" glueballs rather than with 
somewhat artificial static sources. It is a direct application of the Lagrange mesh procedure 
that we presented in Ref. [22]. 
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Our paper is organized as follows. In sec. [TTl we recall the main lattice QCD results 
concerning the lightest scalar glueball. Then, we describe the method to compute the 
effective gluon-gluon potential in sec. IHIl and we comment our results in sec. IIVI We finally 
draw some conclusions in sec. |VJ 



II. RESULTS FROM LATTICE QCD 



An SU(3) lattice calculation in glueball spectroscopy shows that the lightest glueball is 
a scalar particle, whose quantum numbers are J PC = ++ , and whose mass is given by 
1.710 ± 0.130 GeV [l9|]. Theoretical arguments also support this point |25j. The SU(2) 
wave function of this scalar glueball has been first computed in Ref. 20j, and its mass was 
found to be around 1.2 GeV, which is lower than the currently accepted SU(3) value. More 
recently, the SU(3) scalar glueball wave function has been computed [21J] . In this last work, 
it is found that 

m ++ = 1.680 ± 0.046 GeV, (1) 



in agreement with the result of Ref. 



19] . 



The ++ radial wave function which is computed in Ref. 2l| seems thus to be a reliable 
result. As only a few points of this wave function are available, it is more convenient for 
latter calculations to fit them by the trial function 

B~ 



R(r) = exp 



-A 



r 



(2) 



The size parameter r = 0.29 fm = 1.472 GeV -1 is interpreted as the glueball radius in 
Ref. jlll]. Its introduction allows to deal with dimensionless fit parameters A and B. A fit 
with the Levenberg-Marquardt algorithm gives 



A = 0.883 ± 0.045, B = 1.028 ± 0.132, 



(3) 



with a satisfactory agreement since the coefficient of determination is equal to 0.958, close 
to the optimal value of 1. The result is plotted in Fig. [TJ Let us note that R(r) is normalized 
in such a way that -R(O) = 1. 
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FIG. 1: Plot of the ++ radial wave function, taken from Ref. 2l| (circles), with R(0) = 1 and 
ro = 1.472 GeV -1 . The fitted radial wave function ([2]) is also plotted (solid line). 



Consequently, the lattice QCD calculations of Ref. [2l| provide us not only with the ++ 
glueball mass m ++ (PQ), but also with a radial wave function, denoted as R(r) and given by 
Eqs. © and ©. 



III. EFFECTIVE POTENTIAL 



A. Main equation 

In this work, we want to find a Hamiltonian H such that 

HR(r) = m 0++ R(r), 



(4) 



where R(r) and m ++ come from lattice QCD. This eigenequation is the main equation to 
be solved in the following. To our knowledge, this problem has never been addressed before. 
It is of particular interest in order to check the relevance of the approaches with constituent 
gluons. 

In terms of potential models, a glueball can be a bound state of two or more constituent 
gluons coupled in a color singlet. We point out that all the potential models agrees with the 

nnnnn 

fact that the lowest ++ should be a mainly two-gluon bound state |4j, |5|, |6|, |7|, [8|. Moreover, 



the quasiparticle picture of Ref. |3j, based on the QCD Hamiltonian in Coulomb gauge, also 
favors a bound state of two constituent gluons for the lowest ++ . We will thus assume for 
H a two-body standard form, i.e. 

H = T(p 2 ) + V(r), (5a) 

with a local potential V(r) and a semirelativistic kinetic term given by 

T(p 2 ) = 2 ^p 2 + ml (5b) 

Once the mass of the gluon, m g , as well as the radial and orbital quantum numbers n and 
I are specified in = I = for the lightest scalar glueball), the only unknown quantity in 
Eq. (jlj) is the local central potential V(r), that we will compute in the following. 



B. Interpretation of the wave function 

It is worth underlining an important assumption that we make in this work: Motivated by 
the standard representation of the ++ glueball as a two-gluon bound state, we identify the 
lattice wave function R(r), computed from plaquette operators, with the gluon-gluon com- 
ponent of the scalar glueball wave function. This hypothesis can be intuitively understood 
as follows. 

A completely relativistic wave equation for a two-body system is provided by the well- 



known Bethe-Salpeter equation 23]. In this formalism, the relativistic wave function x(r) 
(also known as the Bethe-Salpeter wave function) for a bound state \ fl) of two particles a 
and (3 is defined as 2j| 

x (r) = e- jPR (0| a Q (x a ) a p (x a + r) \Q) . (6) 

In this last equation, |0) is the vacuum state and a 7 is the creation operator for a particle 
7. Moreover, x 7 is the spacetime coordinate of particle 7, r = x@ — x a , and P and R are the 
4-momentum and the spacetime coordinate of the center of mass of the system respectively. 
The relativistic wave function ([6]) can easily be linked to a Schrodinger-like wave function. 
Let us indeed work in the rest frame of the bound state Then, P = (M^,0). In the 
usual nonretarded approximation x° = x°, one has 

x ( r ) = e- iMnt (0| a a (x a ) a p (x a + r) \9) . (7) 



Since the time coordinate t only appears as a global phase, only the spatial separation r, 
computed in the rest frame of the configuration, is now relevant. 

Formula (JTj) can then straightforwardly be applied to the gluon-gluon component of a 
glueball state \G). Indeed, if denotes the creation operator of a gluon, the corresponding 
Bethe-Salpeter wave function at t = is given by 20j 



X (r) = (0| sT J dr Y £m (r) A\{x) A v {x + r)\G) 



(8) 



where the summation on s^ u and the integration on the angular part r enforce a particular 
spin and angular symmetry of the wave function, following the state which is considered. 
Equation (jSJ) can thus be identified as the gluon-gluon wave function for a particular sta- 
tionary glueball state, in the rest frame of the system. 

Interestingly, the wave function (JSj) is precisely what is computed in the lattice QCD 



studies of Refs. 



m 



211 ] for the lightest glueball states. The key point to perform such a 



calculation is the evaluation of the two-gluon operator A^(x) A v (x + r), which is achieved 
on the lattice by computing the correlation matrix between two plaquettes at the different 
points x and x + r. From this discussion, we can justify the identification of the lattice 



wave function of Ref. 
two gluons. 



2l| with the Schrodinger-like wave function of a bound state made of 



Numerical method 



Several methods exist to compute the equivalent local potential from a given wave func- 
tion and its corresponding energy (see Ref. 26|). But, these methods are only applicable 
to the case of nonrelativistic kinematics. We gave in Ref. 22J a procedure, relying on the 
Lagrange mesh method, to make such computations with a semirelativistic kinematics of the 
form (15bj) . We recall here the main points of this method, but refer the reader to Ref. 22] 
for a detailed study. 

The Lagrange mesh method is a very accurate and simple numerical procedure to com- 
pute eigenvalues and eigenfunctions of a two-body Schrodinger equation [27], [28] . It is also 



applicable to a semirelativistic kinetic operator, i.e. the spinless Salpeter equation |29j. In 
the case of radial equations, a Lagrange mesh is formed of N mesh points Xi which are the 
zeros of the Laguerre polynomial L^(x) of degree N {^J. The Lagrange basis is then given 
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by a set of N regularized Lagrange functions, 

fi{x) = (-1)X" 1/2 x{x - x^Lxix) e~ xt \ (9) 

satisfying the condition fj(xi) = \ Sij [27[ and fi(0) = 0. The weights Aj are linked to 
the mesh points Xi through a Gauss quadrature formula which is used to compute all the 
integrals over the interval [0, oof, that is 

poo N 

/ g(x)dx « S*] X k g(x k ). (10) 
Jo k=i 

The regularized wavefunction, given by u(r) = rR(r), is then developed in the Lagrange 



basis. The semirelativistic kinetic matrix elements for the operator f!5bj) can be accurately 
computed in this basi s |29l| . 

We showed in Ref. 22j that, starting from a given regularized wave function u(r) = r R(r), 
obtained in our case from Eq. (J2J), and its corresponding energy mo++, the equivalent local 
potential is accurately given at the mesh points by 

1 - 

V(hxi) = m ++ — s/Tij\/% u{hxj). (11) 

VAi u(hxi) ^ 

In the above equation, h is a scale parameter chosen to adjust the size of the mesh to the 
domain of physical interest. The angular orbital momentum t has to be a priori specified 



22( | . As we deal with the lightest 



in Eq. since the matrix elements T^- depend on 
glueball, we assume £ = and S = in order to obtain a ++ state. Equation (fTTl) requires 
the knowledge of the wave function at the mesh points hxi. So, it is more efficient to work 
with the trial wave function (j2J) than with interpolated points between the few available 
lattice data. 

The first numerical parameter is the number of mesh points, N, that we set equal to 100 



(although iV = 30 already gives a good picture of the potential |22j). The second parameter 
is the scale parameter h. In Eq. ffTTj) . we use a dimensionless variable x, with r = hx. 
A relevant value of h is obtained thanks to the relation h = r a /xN, where xn is the last 
mesh point and r a is a physical radius located in the asymptotic tail of the wave function. 
This radius has to be a priori estimated, but not with a great accuracy, since the method 



is not variational in h 



22 



291 ] . To determine h, we impose that R{hx^)/ max [i?(r)] = e, 
with R{r) given by Eq. (j2J), and e a small number that we arbitrarily fix at e = 10~ 3 . This 
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way of estimating h has already given good results [22]. As it can be seen from Eq. 
max [R(r)] = 1, and we obtain 



x N 



In e 



l/B 

(12) 



A 

Let us note that, strictly speaking, the effective potential V(r) is not unique. Indeed, if 
|0) is a state such that (0 \R) — with (r \R) = R(r), then the whole class of potentials 
U{r) = V(r) + A |0) (0| satisfies the eigenequation (jl]), A being an arbitrary number. But, 
these potentials are clearly nonlocal. Up to our knowledge, there is no physical indication 
that the interaction between two gluons could be non local. We thus compute the unique 
local potential V(r) for the ++ state that we are studying. Consequently, the only possibility 
is A = 0, for which U(r) = V(r). 

IV. RESULTS 

A. Massless gluons 

We begin by setting m g = in the kinetic operator f!5bl) . This corresponds to a model in 
which the gluons are assumed to be massless. With the scale parameter defined by Eq. (|T2l) . 
we can apply formula ( TTT1) . Results are plotted in Fig. El 



The potential obtained by using the optimal values of rriQ++ , A, and B clearly exhibits 
a confining long-range part, and a rapidly decreasing short-range part. The errors on these 
three parameters (see Eqs. ([1]) and ([3])) allow the "true" potential to be located between 
two extremal curves. We have checked that the curves obtained remain stable for different 
values of e, N, and h. These numerical results are compatible with the following Cornell 
potential 

V c (r) =Car- — + D. (13) 
r 

In this expression, as is the strong coupling constant, that we set equal to 0.20 as in the 



case of the static quark-antiquark potential la, p. 42]. The 3 factor is the color factor 



coming from the one gluon exchange between two gluons when the pair is in a color singlet. 
In the confining part, a is the fundamental quark-antiquark flux tube energy density. It is 
usually assumed to be around 0.19 GeV 2 [la . p. 9]. The constant C indicates the scaling 



1 2 3 4 5 6 

r (GeV -1 ) 



FIG. 2: Plot of the optimal potential given by Eq. (fTTI) (dotted line) for m g = 0. It is computed 
from the optimal wave function u(r), plotted with an arbitrary normalization (dashed-dotted line). 
The errors on the glueball mass, mQ++, and on the wave function parameters A and B actually 
allow every potential which is located in the gray area. These results are compared with two 
Cornell potentials CED for the standard values a = 0.19 GeV 2 and a s = 0.20: C = 9/4, D = 
(solid line); C = 3/2, D = 0.30 GeV (dashed line). 

of the energy densrty which is different for a glnon-glnon sys tenr or a q nark-anti q uark pair. 
It is generally assumed to be given by 9/4 (Casimir scaling) [30], or by 3/2 (square root of 



Casimir scaling) [3JJ. The last constant, D, is used to fit the height of potential (1131) on the 
numerically computed optimal potential. 

The Cornell potential with C = 3/2 is closer to the optimal curve than the one with 
the Casimir scaling. However, the interaction with C = 3/2 demands that D = 0.30 GeV. 
This positive constant is somewhat surprising. Indeed, one expects that the parameter D 
encodes some relativistic corrections to the Cornell potential. Talking of such relativistic 
corrections, like spin-orbit or spin-spin terms, is meaningful even in the case of massless 
gluons, since these corrections are expressed in terms of powers of l//i 2 , with a rather large 
value for /j, = (^a/p 2 ^, the dynamical gluon mass generated by the confinement 3^]. Indeed, 
a computation of /i with the fitted wave function leads to /i — 0.585 ±0.128 GeV. Since 
we consider here a £ = S = state, the only remaining relativistic corrections are retardation 
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terms and contact (spin-spin) interactions, which are both negative in this case [4], |33j . The 
Casimir scaling seems thus more satisfactory since it is compatible with a vanishing value 
of D and it is still located in the allowed region. 

One can see in Fig.[2]that the very short-range part of the numerically computed potential 
is less singular than the Coulomb potential. This part is actually very sensitive to the short- 
range behavior of the wave function. It has been shown that the wave function of a true 

n 

semirelativistic Coulomb problem diverges in r = [34J], but such a divergence cannot be 
computed in lattice calculations. Anyway, more points should be necessary to elucidate this 
short-range behavior. 

Let us note that the ground state masses of the Hamiltonian (J5j) with the two fitted 
potentials ( 1131) are respectively 1.675 GeV and 1.670 GeV with C = 9/4 and C = 3/2. These 
values are contained in the error bars of mo++ . The corresponding wave function with C = 
3/2 is nearly indistinguishable of the optimal wave function ([2D, while the one for C = 9/4 
decreases more quickly. This is coherent with Fig. [2j where we can see that the Cornell 
potential with the Casimir scaling increases faster than the optimal curve, consequently 
leading to a faster decreasing wave function. 

B. Massive gluons 

As mentioned in the introduction, it is argued in many works that the constituent gluons 
should have a fixed nonzero mass, typically around 0.5 ± 0.2 GeV (5), [ljj. The effective 
gluon-gluon potential can also be computed by using our method. In order to clearly see 
the changes between m g = and m g > 0, we have computed the effective potential for a 
rather large but still relevant value of m g , namely 0.7 GeV. The result is plotted in Fig. [3j 



A comparison with Fig. [2] shows that the asymptotic behavior of V(r) drastically changes 
when m g > 0. Instead of a monotonically increasing potential, it now seems to saturate at 
some energy scale. Such a saturation suggests an interpretation in terms of string breaking 
effects, which have been measured in lattice QCD in the case of a static quark- ant iquark 
pair. Thus, the optimal potential of Fig. [3] should no longer be compared to a Cornell form. 
However, in a pioneering work about glueballs as bound states of massive gluons by Cornwall 
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FIG. 3: Plot of the optimal potential given by Eq. ([lip (dotted line) for m g = 0.7 GeV. It is 
computed from the optimal wave function u(r), plotted with an arbitrary normalization (dashed- 
dotted line). The errors on the glueball mass, m,Q++, and on the wave function parameters A and 
B actually allow every potential which is located in the gray area. These results are compared with 
two "Cornwall and Soni potentials" of the form flUD with S 2 = 0, r c = 2m g /Ca, a = 0.19 GeV 2 , 
and a s = 0.67: C = 9/4, D = -0.67 GeV (solid line); C = 3/2, D = -0.50 GeV (dashed line). 

and Soni it is suggested that the gluon-gluon potential at lowest order reads 



where S is the spin of the glueball. This potential clearly saturates at a value 2m g + D } then 
forbidding bound states with a mass greater than 4m g + D. For r r c , the confining part 
reduces to (2m g /r c ) r, and a linearly rising confining potential is recovered. By comparison 
with the Cornell potential (|T3|) . one could expect that 



which provides us with a definition of r c in terms the more intuitive parameters m g and 
a. Some useful remarks can also be made about the short-range part of the potential (1141) . 
i.e the part which is proportional to a s . Firstly, it involves a Yukawa potential instead 
of a Coulomb one, because both exchanged gluons and constituent gluons are massive. 




(14) 
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Moreover, it is spin-dependent even at the lowest order. This is a characteristic feature of 
effective potential with massive gluons. Other spin-dependent corrections in 1/m 2 , also exist, 
but we do not mention them since they are absorbed in the constant D by definition. Again, 
this constant is fitted so that the absolute height of the optimal potential is recovered. 

It is readily observed in Fig.[3]that potential (1141) fits rather well the numerically computed 
potential, in particular the long-range behavior, for r > 2 GeV" 1 . For comparison, we notice 
that r c = 3.27 GeV" 1 when C = 9/4 and 4.91 GeV" 1 when C = 3/2. But, as we argued in 
the previous section, the short-range part is only poorly known with the current lattice data. 
The difference between C = 9/4 and 3/2 is less important than for the Cornell potential even 
if the Casimir scaling seems closer to the optimal curve when the usual value a = 0.19 GeV 2 
is used. It is worth mentioning that we fitted a s = 0.67, which is higher than in the case of 
massless gluons. But, such a value is commonly used in models with massive gluons 

We finally note that the ground state masses of the Hamiltonian §5§ with the two fitted 
potentials (EH]) are respectively 1.724 GeV and 1.721 GeV with C = 9/4 and C = 3/2. Again, 
these values are contained in the error bars of rriQ++. In this case, the corresponding wave 
functions with both values of C are very similar to each other, but with a slightly larger 
spatial extension than the optimal wave function ([2]). 

C. Mass dependence of the potential 

A plot of the effective potential for different values of m g is given in Fig. 0] in order to 
see more clearly the evolution of the potential shape with an increasing gluon mass. It is 
readily observed that there is a transition between a linearly rising regime at zero mass and 
a saturating one at high gluon mass. 



As we already pointed out, the short-range behavior of the computed potentials is rather 
poorly known because of the current precision of the lattice data. However, interesting 
informations can be deduced concerning the long-range part of the potential, typically for 
r > 3 GeV -1 . One can find indeed that the different potentials computed in Fig. H]can be 
correctly described at large r by the following form 

V f (r,m g ) = — (l-e-^ m * r )+D, (16) 
7m 9 ' 
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FIG. 4: Plot of the optimal potential given by Eq. for different values of m g (black lines). 
Potentials are computed from the optimal values of the wave function u(r) and of the mass m ++ • 
The long-range potential Vy(r, m g ), given by Eq. (|16p. is plotted for the different values of m g (gray 
lines) with C = 3/2, a = 0.19 GeV 2 , and 7 = 0.362. The constant D is fitted so that the absolute 
height of the corresponding numerically computed potential is well reproduced. 

with 7 = 0.362 ± 0.021 a fitted parameter. Such a form has the advantage of having a 
non trivial limit for m g = 0, which is simply the linearly rising potential Car + D. When 
m g > however, it saturates at the value (Ca/-ym g ) + D. Potential ffTBl) has been plotted 
for the different values of m g that we considered, with the previously found values of <r, 
that is 0.19 GeV 2 . We fixed C = 3/2 because we already saw in sec. IIV Al that the optimal 
numerically computed potential was better reproduced with that value of C. But, again, we 
stress that the precision on the potential does not allow to decide whose scaling law is the 
best one. The constant D have been fitted for each m g in order for the potential ( TL6l) to 
have the right absolute height. 

In the previous section, we showed that the potential obtained with m g = 0.7 GeV can 
be well reproduced by the form (jl4| . whose long-range part is a priori inequivalent to the 
potential ( |T6l) . Actually, both expressions coincide if 




(17) 
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With a = 0.19 GeV 2 , this relation states that potentials (1141) and (1161) are identical for 
m g = 0.768 ± 0.024 GeV when C = 3/2, and m g = 0.627 ± 0.020 GeV for C = 9/4. 
Remarkably, this corresponds to a gluon mass around 0.7 GeV, that is the case that we 
treated previously. 

V. CONCLUSIONS 

In this work, we have shown how to compute the effective potential between two gluons 
from the mass and wave function of the ++ state obtained in lattice QCD, reasonably 
assuming that the lattice wave function mainly gives the gluon-gluon part of the glueball 
Schrodinger-like wave function. This method, which is here used for the first time, has the 
advantage of allowing to deal with semirelativistic kinematics, which is necessary to include 
massless gluons systems in the discussion. 

The relevance of potential models to describe gluons is still a matter of controversy 
nowadays. In particular, the mass of a constituent gluon is an open problem: Should it 
be zero or not? Arguments favoring both hypothesis can be found in the literature, so we 
computed the effective gluon-gluon potential in both cases. When the constituent gluons 
are massless, the potential we find is compatible with a Cornell one for standard values of 
the parameters. When the constituent gluons are massive however, the effective potential is 
merely compatible with a saturating one with a Yukawa- type short-range part. 

These results go beyond the usual computation of the static potential between two static 
sources in the adjoint representation. In this work indeed, we started from pure gauge results 
in lattice QCD and we showed that the constituent gluon picture naturally emerges from 
these lattice data. This could yield an a posteriori justification of the success of potential 
models in the description of the glueball spectrum. Moreover, we showed that the allowed 
shape of the potential is constrained by the mass of the constituent gluon (see Eq. ( fl6l) ): 
The linear confinement or the saturating one are only valid for massless or massive gluons 
respectively. Although the Casimir scaling seems favored, the current precision of the lattice 
data does not allow to decide which scaling law is correct. This problem could be solved 
by making high precision computations of the wave function. We mention however the 
recent work of Ref . [35J , where the behavior of a quark-antiquark-constituent gluon system 
is investigated on the lattice, with results in agreement with the Casimir scaling hypothesis. 
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An important application of our method is that, provided accurate wave functions for 
glueballs with a higher total spin are computed (too few points are currently available 



for the 2 ++ wave function 21|). spin-dependent terms-vanishing for the scalar glueball- 
could be studied. This would be a check of the relevance of approaches with massless or 
massive constituent gluons. Indeed, potentials coming from both models already differ at 
the lowest order for what concerns spin-dependent terms. We consequently think that a 
very accurate determination of the wave functions of the lightest ++ and 2 ++ glueballs 
in lattice QCD, especially at short range, could serve as a check to determine the most 
relevant potential approach of glueballs. We expect that, provided these wave functions are 
known, our method should lead to numerically computed potentials with a common long- 
range behavior-confinement is indeed not supposed to depend on J PC . On the contrary, 
the short-range behavior should be different, and give important informations about the 
relativistic corrections. We hope to present such a study in future works. 

Acknowledgments 

The authors thank the FNRS for financial support, and Francis Michel for stimulating 
discussions. 



[1] H. B. Meyer, |hep-lat/050 8002 ; W. Ochs, |hep-ph/0609207| 

[2] C. E. Carlson, T. H. Hansson, and C. Peterson, Phys. Rev. D 27, 1556 (1983). 

[3] A. Szczepaniak, E. S. Swanson, C.-R. Ji, and S. R. Cotanch, Phys. Rev. Lett. 76, 2011 



(1996) [hep-ph/9511422 ; A. Szczepaniak and E. S. Swanson, Phys. Lett. B 577, 61 (2003) 



|>ep-ph/0308268]. 

[4] H. Fritzsch and P. Minkowski, Nuovo Cimento 30 A, 393 (1975); T. Barnes, Z. Phys. C 10, 
275 (1981). 

[5] J. M. Cornwall and A. Soni, Phys. Lett. B 120, 431 (1983). 

[6] A. B. Kaidalov and Yu. A. Simonov, Phys. Lett. B 636, 101 (2006) [hep-ph/051215l] . 
[7] F. Buisseret, arXiv:0705.0916. 

[8] F. Brau and C. Semay, Phys. Rev. D 70, 014017 (2004) |hep-ph/04 12173] , and references 



16 



therein; V. Mathieu, C. Semay, and B. Silvestre-Brac, Phys. Rev. D 74, 054002 (2006) 
[hep-ph/0605205] . 

F. Buisseret and C. Semay, Phys. Rev. D 73, 114011 (2006) [hep-ph/0511270] . 

F. J. Llanes-Estrada, P. Bicudo, and S. R. Cotanch, Phys. Rev. Lett. 96, 081601 (2006) 

|hep-ph/0507205] . 



F. J. Llanes-Estrada and S. R. Cotanch, Phys. Rev. Lett. 84, 1102 (2000) [hep-p h/9906359|. 
W.-S. Hou and A. Soni, Phys. Rev. D 29, 101 (1984); W.-S. Hou, C.-S. Luo, and G.-G. Wong, 
Phys. Rev. D 64, 014028 (2001) |hep-ph/010iT46] ; F. Brau and C. Semay, Phys. Rev. D 72, 



078501 (2005) [hep-ph/0411 108|. 

A. C. Aguilar, A. Mihara, and A. A. Natale, Phys. Rev. D 65, 054011 (2002) [hep-ph/0109223] . 
N. Brambilla, P. Consoli, and G. M. Prosperi, Phys. Rev. D 50, 5878 (1994). 
G. S. Bali, Phys. Rep. 343, 1 (2001) |hep-ph/0001312|. 



Yu. A. Simonov, |hep-ph /9911237 A. B. Kaidalov and Yu. A. Simonov, Phys. Lett. B 477, 
163 (2000) [hep-ph/9912434] . 

G. S. Bali et al, Phys. Rev. D 71, 114513 (2005) [hep-lat/0505012] . 

J. Vijande, P. Gonzalez, H. Garcilazo, and A. Valcarce, Phys. Rev. D 69, 074019 (2004) 
hep-ph/0312165|; P. Gonzalez, A. Valcarce, J. Vijande, and H. Garcilazo, Int. J. Mod. Phys. 
A 20, 1842 (2005) [hep^ph/0409202 . 



C. J. Morningstar and M. Peardon, Phys. Rev. D 60, 034509 (1999) |hep-lat/9 901004|; Y. 

Chen et al, Phys. Rev. D 73, 014516 (2006) [hep-lat/0510"074] . 

P. de Forcrand and K.-F. Liu, Phys. Rev. Lett. 69, 245 (1992). 

M. Loan and Y. Ying, Prog. Theor. Phys. 116, 169 (2006) [hep-lat/0603030 1. 



F. Buisseret and C. Semay, Phys. Rev. E 75, 026705 (2007) [physics/0611146]. 

E. E. Salpeter and H. A. Bethe, Phys. Rev. 84, 1232 (1951); J. Schwinger, Proc. Natl. Acad. 
Sci. 37, 455 (1951). 

G. C. Wick, Phys. Rev. 96, 1124 (1954); R. Cutkosky, Phys. Rev. 96, 1135 (1954). 
G. B. West, Phys. Rev. Lett. 77, 2622 (1996) i ;hep-ph/9603"31~6] . 

A. Lovell and K. Amos, Phys. Rev. C 62, 064614 (2000) [nucl-th/0007031 1, and references 
therein. 

D. Baye and P.-H. Heenen, J. Phys. A 19, 2041 (1986). 

M. Vincke, L. Malegat, and D. Baye, J. Phys. B 26, 811 (1993); D. Baye, J. Phys. B 28, 4399 



17 



(1995). 

[29] C. Semay, D. Baye, M. Hesse, and B. Silvestre-Brac, Phys. Rev. E 64, 016703 (2001). 



[30] G. S. Bali, Phys. Rev. D 62, 114503 (2000) [hep-lat/0006022 1; S. Kratochvila and P. de 



Forcrand, Nucl. Phys. B 671, 103 (2003) [hepd at/0306011 1; C. Semay, Eur. Phys. J. A 22, 

353 (2004) |hep-ph/0409"lQ5] . 
[31] K. Johnson and C. B. Thorn, Phys. Rev. D 13, 1934 (1976); G. Martens, C. Greiner, S. 

Leupold, and U. Mosel, Phys. Rev. D 70, 116010 (2004). 
[32] A. B. Kaidalov, Yu. A. Simonov, Phys. Lett. B 477, 163 (2000) [hep-ph/ 991 2434] . 



[33] V. Mathieu and F. Buisseret, hep-ph/0702226 



[34] B. Durand and L. Durand, Phys. Rev. D 28, 396 (1983); 50, 6662(E) (1994); P. Castorina, P. 

Cea, G. Nardulli, and G. Paiano, Phys. Rev. D 29, 2660 (1984). 
[35] P. Bicudo, M. Cardoso, and O. Oliveira, hep-ph/0704.2156. 



18 



